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Abstract — While nonlinear optical spectroscopy is becoming 
more commonly used to study the excited states of nonlinear- 
optical systems, a general theory of inhomogeneous broadening 
is rarely applied in lieu of either a simple Lorentzian or Gaussian 
model. In this work, we generalize all the important linear and 
second-order nonlinear susceptibility expressions obtained with 
sum-over state quantum calculations to include Gaussian and 
stretched Gaussian distributions of Lorentzians. We show that 
using the correct model to analyze experiments that probe a 
limited wavelength range can be critical and that this theory is 
better able to fit the subtle spectral features - such as the shoulder 
region of a resonance - when both models produce qualitatively 
similar responses. 



I. Introduction 

WHEN calculating the nonlinear-optical susceptibilities 
near resonance, it is important to account for damping 
to get a reasonable approximation to the dispersion of a real 
molecule. Clearly, neglect of damping is catastrophic since 
sum-over-states (SOS) theories predict an infinite susceptibil- 
ity on resonance. The SOS expression, as derived by Orr and 
Ward,[l] naturally admit damping corrections by making the 
eigenenergies of the system complex. The imaginary part of 
the energy is related to the resonance width and is called 
the natural linewidth. A stationary isolated molecule's line- 
shape is thus described by a Lorentzian. When the width is 
related to the intrinsic properties of a molecule, the width of 
a peak for a collection of noninteracting molecules is referred 
to as homogeneous broadening, i.e. the observed spectrum is 
the same as that of a single molecules times the number of 
molecules. 

If the molecules interact or are placed in an environment 
with random perturbations, the natural linewidth is smeared 
out, leading to what is called inhomogeneous broadening 
(IB). Stoneham described in a detailed review article how the 
linear absorption spectrum in a crystal could be accurately 
modeled by applying stochastic averaging to a distribution 
of Lorentzians. |2| Since such statistical processes often lead 
to gaussian peak shapes such as the Boltzmann velocity 
distribution, systems with resonance features whose widths are 
much larger than the natural line width are often modeled by 
Gaussian functions. 
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Dirk et al[3 | and Berkovic et al[4| showed that it was impor- 
tant to take damping into account by using a complex energy 
in the Orr and Ward SOS expression even when considering 
off-resonance susceptibilities. The distribution of molecules, 
when embedded in a polymer, are randomly perturbed by 
virtue of the fact that each molecule experiences a different 
local environment due to microscopic inhomogeneities. The 
distribution of sites around a molecule in a polymer can be 
measuredj5), and quantified by a stretched exponential 
distribution function. 

Many excited states are often needed to correctly de- 
scribe the dispersion of the hyperpolarizability [7|, [8| and 
vibronic overtones may complicate the analysis. [9], [10| Since 
it has been shown that vibronic states are unimportant off- 
resonance, |flj], 132] and on-resonance, we assume that their 
affect is smeared out due to thermal broadening (i.e. thermal 
fluctuations and the distribution of sites cause variations that 
are larger than vibronic energies), we argue that neglect 
of vibronic states and applying a stochastic averaging over 
Lorentzian functions using a stretched exponential weighting 
distribution function leads to a good model of the nonlinear 
spectra. 

In this work, we focus on using the SOS expression to 
calculate the linear and second-order nonlinear-optical re- 
sponse under inhomogeneous broadening using the above 
approximations. (We will apply the same approach to the 
dipole-free SOS expression] 13 1 in a future publication.) Our 
theoretical results are compared with experimental data in dye 
doped polymers for validation. 

II. Theory 

Sum-over states quantum perturbation treatments of the 
b th -oxdex nonlinear susceptibility tensor, £^ fc , in the dipole 
approximation yields a sum of terms of the form:[T| 
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where (/x») nm is the nm-matrix element of the i th component 
of the electric dipole operator, 0J nm the transition frequency 
(energy) between states n and m, u>i the frequency of the i th 
optical field, and r„ g the phenomenological damping factor. 
The numerator is a product of b + 1 transition moments and 
the denominator a product of b energy terms. For an isolated 
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molecule, the damping factor T ng is inversely proportional to 
the lifetime of state n and is a measure of the width of the peak 
in the spectrum of Qj k associated with a transition between 
state n and the ground state g. 

In real systems, molecules interact with each other yielding 
a broadening of the peaks in a spectrum. One common method 
to treat this case is to adjust the parameters T ng , If the statistics 
of the molecules are Gaussian (as one finds in Doppler 
broadening), the Lorentzian function resulting from adjusting 
T ng has the correct width but does not have a Gaussian shape. 
Another common method for modeling a nonlinear spectrum 
is to assume the peaks are Gaussian in shape. Depending on 
the region of interest (i.e near the peak or in the tail of the 
spectrum) one model may be more applicable than the other. 

An exact method of treating inhomogeneous broadening 
is to apply the statistics of the molecular interactions to the 
susceptibility. Each molecule in an ensemble is then viewed 



as having a different transition frequency (energy), 
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where Suj ng — uo ng — uj ng , uj ng is the mean value of the 
transition frequency, Nfjjng) the normalization factor, and 7„ g 
the linewidth of the distribution - the susceptibility will be of 
the form, 
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Toussaere developed such a 
generation (SHG) and the Pockels effect using Gaussian 
statistics.lfl4l 

In materials such as dye-doped polymers, the statistics that 
best model the effect of the distribution of sites on processes 
such as relaxation of molecular orientation order (which is 
normally an exponential process as described by Debye) is 
a stretched exponential. Using this analogy, we propose that 
variations in the local electric fields in a polymer yields 
the same statistics for modeling the susceptibility, so the 
distribution is also assumed to be a stretched Gaussian of the 
form, 
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where (3 is the distribution of sites parameter, and for many 
systems - such as in dye-doped polymers - varies from (3 = 
for an infinitely broad distribution to (3 = 1 for a single 
characteristic width. Again, AT (7, 0) is a normalization factor, 
which depends also on (3, and will be written as 

iV(7„ 9 ,/3)=7n 9 V^(/3), (5) 

where, 
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to remain compatible with previous inhomogeneous broaden- 
ing representations which use Gaussian statistics fl4l . |fl5ll , 
|fl6| , |[T7ll . Indeed, such statistics other than pure Gaussian 
were observed in light scattering experiments ifTHl . |[T6ll . ifTTll . 

A. Normalization and Approximations 
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Fig. 1. Relative difference between B(f3) when the integral is taken over 
all frequency space, and when the lower limit for B(/3) is given by ui ng = 2 
eV 

Because B(f3) is difficult to calculate, we need to make 
simplifying assumptions. Typical inhomogeneous linewidths 
are narrow compared with the transition frequency |fl9l , IfTTll . 
so we can integrate over all space even though there is a finite 
lower limit to the integral. Furthermore, we assume that B(f3) 
is independent of 7„ 9 . Figure Q] shows the relative difference 
between B((3) determined by integrating of over all frequency 
space and integration to a lower limit of oj ng = 2 eV. Clearly, if 
the inhomogeneous linewidth is smaller than 70 meV - which 
is the case for all our experimental data - this assumption 
leads to a negligible relative difference (< 1%) for (3 > 0.3. 



25 



x 10 



20 -; 



CO. 

5 



CO. 

5 



CO. 

5 



= 30 meV 
= 40 meV i 



Y =60meV g ,„: 
y ng =7()mev| 10 




0.8 0.9 



Fig. 2. Relative difference between B(/3, "(ng) and B({3, 7„ 9 = 50 meV) 



^ The second assumption is that the normalization function 
B{(3) is independent of the inhomogeneous linewidth, 7„ 9 . 
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Figure [2] shows that this assumption is valid to good approxi- 
mation (within 1%) for /3 > 0.15 (see inset in Figure |2). The 
reference value of B((3, j ng = 50 meV) was chosen because it 
is a typical linewidth for the chromophores we consider. These 
assumptions allow us to calculate one normalization function 
independent of u> ng and -f ng when we consider (3 > 0.3, which 
is critical for decreasing the computation time of the theory 
when applied to a least-squares fit to the data. 

In this paper, we derive the expressions for the most 
important linear and second-order susceptibilities for any 
general statistics. Specific results are presented for Gaussian 
and stretched exponential suceptibilities that model the SHG 
response (hyperpolarizability) for a two-level system and data 
from linear absorption experiments (polarizability) on sil- 
icon phthalocyanine-(mono)methylmethacrylate(SiPc) doped 
(poly)methylmethacyrlate (PMMA) thin films. All the expres- 
sions for the second-order processes are displayed in detail in 
Appendix B. 

B. First-Order Susceptibility 

In this section, we outline the approach of the calculation for 
a general distribution function for a linear susceptibility (which 
has one energy denominator) and in the next section apply the 
theory to a second-order susceptibility (which has two energy 
denominators). These two cases can be easily generalized to 
any-order susceptibility and any statistics. 

From time dependent perturbation theory, the 1-D first-order 
molecular susceptibility is lEOl . |[T9l : 

(-w; u) = -AllYl {^ngD^i-w; uj)} , (7) 
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The inhomogeneous broadened Lorentzian function is then 

Di B (-u>;u>)= / D^(-uj;uj)f ng (6uj ng )d(6uj ng ). (9) 



Once £)„ (— is known it can be substituted into Eq. 
(0 for D^(~uj;uj) to create an inhomogeneously broadened 
microscopic susceptibility. To this end, we proceed to deter- 
mine the first term in Equation (O. Substituting Equation (0]i 
for the Gaussian function, changing the integration variable to 
t = (uj' ng — uj ng ) /j n g, and rearranging the denominator so that 
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To simplify the expression in Equation ilOi , we define 

Wg k \z) as the following: 
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and the complex error function as [21], 

exp(-t 2 ) 



W{z) 



t 



dt 



(12) 



The second term in Eq. (O is derived analogously so we 
write the first-order energy denominator that accounts for 
inhomogeneous broadening as: 
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In the limit that = 1, the first-order energy denominator can 
be written in terms of the complex error function, W(z) 1 14|, 

ma, nsi: 
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which can be evaluated for all values of z using the results 
in Abramowitz and Stegun [21 1. The general transformation 
from the Lorentzian model to the inhomogenous broadening 
model is summarized in Table [HI] 

C. Second-Order 

Using time dependent perturbation theory, we define a 
second-order Lorentzian energy denominator as: 

Df m (-u) a ;oJi,u)2) = Si, 2 j[(fiz a - to a )(Cl mg - lui)]^ 1 
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- [(f2* fl +UJ 2 ){Q m g - Wl)] 
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where Si ^ is the permutation operator, which averages over 
all distinct permutations of and u>2- 

It is significantly more difficult to transform nonlinear 
Lorentzian energy denominators to nonlinear IB energy de- 
nominators because of the product of Lorentzian terms in 
the denominator. In order to transform Df m for a specific 
experiment (i.e. for specific input and output frequencies), the 
number of excited states must be known prior to performing 
a partial fraction expansion on each term in the energy 
denominators. For example if there is only one distinct excited 
state (/), it maybe necessary to perform the following partial 
fraction expansion, 
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to eliminate the product of the two Vli g terms. These type of 
expansions allow us to write the nonlinear energy denomina- 
(20) tors i n terms of Wp{z) or complex error functions. 



However, a perfect square in the denominator requires a 
different approach. The approach for evaluating the funda- 
mental transformations from the homogeneous formulation to 
the inhomogeneous formulation for the quadratic dependencies 
(cubic dependencies are studied in Ref. [22|) on the transition 
frequency, uji g is shown below. 
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Beginning with a quadratically dependent term like the 
following, 

C 2 
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we integrate its product with the distribution function (Equa- 
tion ©), 
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as the initial step in the transform. Substituting Equation (|4]i 
for fng(w' ng — io ng ), and changing the integration variable to 
t = (ui' ng — <jJng)/lng, results in the following: 
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where z = (— u> ng + iT ng + uj)/j ng . In the limit that (3 = 1, 
Equation ( fT9l ) looks very similar to Equation ( fTOb just before 
substituting for W(z) except that the denominator in the 
integral is second-order in (z — t). To reduce the denominator 
to first-order in (z — t) so that the integral can be replaced 
with W(z), we perform integration by parts twice: 
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where we have used T = exp(— i 2 ) to simplify the presenta- 
tion. The first term on the right on the top line of Equation (|20| > 
vanishes because the argument of the exponential is w — 10 3 
at the lower limit ( usually uj 7 in the visible). Note that 
Equation d20b is kept general (the lower limit of integration is 
not approximated as —00). We only use this approximation as 
needed for special cases. To get the second line of Equation 
(|20| |. we use t/(z — t) = z/(z — t) — l. Therefore we can write 



Equation ( TT8l . when (3 = 1, as the following, 
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Thus we have derived the convolution of a second-order 
denominator term with a Gaussian, which can be generalized 
to any complex term that has a second-order dependence on 

(z-t). 

Tableimiin Appendix A summarizes the fundamental energy 
denominators up to second-order for the Lorentzian and IB 
theories with (3 < 1, and (3 — 1, respectively. These two 
transformations can be used to construct any IB energy de- 
nominator for any first-, and/or second-order response. 

TABLE I 

Compact form ovWi x '(z) and W(z). 
(3 Term Compact Form 
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The energy denominators Z?; m ... for the secondhand 
higher)-order susceptibilities are complex combinations of 
Wg X \z) or W(z). To make them more readable we have 
developed a compact notation. For j3 < 1, we have added 
a subscript to (3 to indicate the excited state involved in the 
process and a * on the power of W to indicate a complex con- 
jugate of the complex argument fl ng = uj ng — iT ng . Similarly 
for (3—1, the subscript on W indicates the excited state and 
the superscript * on W indicates the complex conjugate of 
fl This allows us to use a simple frequency argument of ±uj 
which significantly improves the readability of the equations in 
Appendix B. An example of the transformation to the compact 
notation is given below: 
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This compact form is sufficient to describe all of the inhomo- 
geneous broadening contributions to the energy denominators 
in this paper because the arguments are all of the form 



or an( j a u f t ne excited state transitions 

are from/to the ground state (g). Table |I| summarizes the 
transformations to the compact notation. 



D. Linear Susceptibility 

In this section we discuss the effects of (3 on the linear 
response of a dye-doped polymer system with one excited 
state. Figure [3] shows the real and imaginary parts of the linear 
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Fig. 3. Real and imaginary [inset] parts of p^W (— u>; uj) for a single excited 
state centered at 653 nm for various values of j3. For all values of /3, V = 10 
meV and 7 = 50 meV. The homogeneous-broadening (Lorentzian) theory of 
the electronic transition is denoted by L and T = 50 meV. 

susceptibility, which is proportional to Eq. (|7]i with n = 1 and 
the appropriate energy denominators, for various values of (3 in 
the IB model. For comparison the homogeneous (Lorentzian) 
theory is denoted by L. We use typical values for the transition 
moment and homogeneous and inhomogeneous linewidths. In 
this example the excited state is centered at about 653 nm 
with a transition moment of 11.5 D. For all the IB curves 
the Lorentzian linewidth is 10 meV, and the IB linewidth 
is 50 meV. The Lorentzian curve is generated using a 50 
meV linewidth and the same center frequency and transition 
moment as the IB curves. As we expect the magnitude of the 
response for a fixed inhomogeneous linewidth decreases and 
becomes broader as the distribution of sites within the polymer 
matrix increases (decreasing (3). To better describe the effects 
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°f X^H~ as a function of (3 (distribution of sites) and 
7„ g (IB) for fixed r raff (L) and x^H - w;ll>) as a function of 
(3 and T ng for fixed 7„ g . In Figure [4] the single excited 
state transition is centered at 653 nm, with a transition dipole 
moment of 11.5 D, and a fixed Lorentzian linewidth of 10 
meV, while Figure [5] uses a fixed IB linewidth, j ng of 50 
meV. 
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Fig. 5. Contours of the real and imaginary parts of X [~ generated 
using the generalized IB theory for a single excited state centered at 653 nm 
(fJ-ng = 11.5 D, 7n ff = 50 meV) for three probe wavelengths (615 nm, 630 
nm, 645 nm). 

Figures @] and [5] each show a distinct character of 
X^(— w,oS) for wavelengths between 615 nm and 645 nm. 
The (3 parameter therefore gives an additional degree of 
freedom for modeling guest/host materials where a distribution 
of sites requires using non-Gaussian statistics for both the 
index of refraction and the absorption coefficient. 

III. Second-order Susceptibility 

In this section, we demonstrate that there are significant 
differences in the response from two-level systems depending 
on the broadening model. We concentrate on second harmonic 
generation (SHG) as an example because of its near-universal 
use in characterizing molecules. For a two-level system, the 
SHG hyperpolarizability ( — 2w; ui, uj) oc Da{— 2ui; u, ui) 
[see Eq ( 1231 and Eq. d25l l for the IB models, and Eq. ( fT~5T > for 
the L model] where I in this case is 1. 

TABLE II 

Inhomogeneous and Lorentzian broadening parameters for a 
two-level system. 



Model wig Tig 7; fl 

L 1.85 eV 0.05 eV 

IB 1.85 eV 0.025 eV 0.05 eV 



Fig. 4. Contours of the real (left) and imaginary (right) parts of x' 1 ' (— ui; ui) 
generated using the generalized IB theory for a single excited state centered 
at 653 nm (p n g = 11-5 D, T ng = 10 meV) for three probe wavelengths 
(615 nm, 630 nm, 645 nm). 

of the distribution of sites, IB, and homogenous broadening 
in the linear response we plot the real and imaginary parts 



Figure [6] shows the response of a two-level system using 
the parameters in Table [TT] For the IB model four values of (3 
have been studied as shown in the legend. It is clear that the 
response in the IR for all models are qualitatively similar but 
quantitatively different. Thus a fit to an SHG experiment that 
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Fig. 6. Second harmonic generation response for several inhomogeneously 
broadened two level systems. The real part of x' 2 ' (—2a;; ui, ui) is shown as a 
function of the incident wavelength. A homogenously broadened (L) response 
is included for comparison. The inset shows the magnitude on a log scale in 
the visible region. Parameters for the curves can be found in Table ITT] 

probes this system in the IR region would come to a similar 
conclusion for the parameters describing which excited states 
are involved in the process. The difference between the two 
models could be used to estimate the error. 

In the visible part of the spectrum (see inset), the different 
models do not agree with each other, where the L model 
is dramatically different both qualitatively and quantitatively 
from the IB models. This becomes important for those exper- 
iments that probe the system over a limited wavelength range 
(or more egregiously, at a single wavelength) and attempt to 
determine an off-resonance value for the hyperpolarizability. 
In this example, an SHG experiment that probed this one- 
level system only in the visible would determine a completely 
different set of parameters depending on the model used to 
interpret the results. 
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Fig. 7. Linear absorption of SiPc/PMMA in comparison to least-squares 
fits using Lorentzian and IB theories ( L={/xi 9 = 8 D, and T\ g = 18 meV 
}, IB(/3 = 1) = {m g = 1A D, F lg = 9 meV, and ji g = 19 meV }, 
IB(/3 = 0.7) = {fiig = 7.5 D, T lg = 7 meV, and -y lg = 18 meV }). 

the low-energy side. This gives a least squares fit of j3 = 1. In 
either case, it appears that the distribution of sites is narrow in 
the SiPc/PMMA system. The result is reasonable given that the 
SiPc molecule is non-dipolar so there are no dipole interactions 
with the surroundings. All broadening therefore results from 
higher-order multipoles, which contribute more weakly to the 
width. 

Given that the calculated homogeneous broadened 
linewidths for each model compare well with the literature 
on the temperature dependence of these linewidths for 
doped PMMA systems,] 24 J it is clear that the IB model 
describes the response better than the Lorentzian model 
especially in the wings of the resonance. It should therefore 
give a better estimate of zero-frequency susceptibilities 
for comparison with the Thomas-Kuhn sum-rule quantum 
limit. [25 1, l26l . ll27l . ll28l Indeed, it would be interesting to 
apply inhomogeneous broadening models to the fundamental 
limits of the dispersion of the hyperpolarizability. ||29l It 
may very well be the case that the gap between the best 
nonlinear-optical molecules and the fundamental limit ifTTI . 
[12], [ 28 1 may be partially explained by inaccuracies in the 
dispersion model. 



IV. Comparison of Theory to Experimental 
Results 

Next, we compare IB Lorentzian models with linear absorp- 
tion data from (SiPc/PMMA) thin films|23| to show that the 
wing region of the linear absorption resonance is often better 
described using an inhomogeneous broadening model. A sim- 
ilar comparison for the quadratic electroabsorption spectrum 
of this molecule can be found in the literature. [23), |22| 

Figure [7] compares SiPc experimental linear absorption data 
with the Lorentzian (L) theory, and the IB theory for two 
values of (3. A least squares fit gives (3 = 0.9. It is possible that 
the high-energy side of the peak is affected by the next excited 
state, so a least squares fit was performed from the peak to 



V. Conclusion 

In conclusion, we calculate the inhomogeneously broad- 
ened linear and second-order nonlinear susceptibilities for a 
Gaussian and stretched Gaussian distribution of Lorentzians. 
A Lorentzian model is found to be inaccurate in predicting 
the shape of the linear absorption spectrum of SiPc/PMMA. 
However, the IB-broadened theory fits the data over a broad 
range of wavelengths and shows that the distribution of sites 
is nearly Gaussian, implying that interactions between the 
polymer and dopant are small, as we would expect of a non- 
dipolar molecule. 

Since the nature of the broadening mechanisms, especially 
for nonlinear spectroscopy, is shown to affect the dispersion 
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dramatically, the determination of excited state properties 
of molecules from spectroscopy requires that such a theory 
be used. Indeed, the inconsistency of transition moments as 
determined by different processes (i.e. linear versus nonlinear 
spectroscopy) may be due in part to the use of inappropriate 
dispersion models. Furthermore, the practice of extrapolating 
single-wavelength measurements to get the off-resonance hy- 
perpolarizability, /3n,ll3CT may lead to large inaccuracies. As 
such, IB theory may be an important tool for interpreting any 
nonlinear-optical experiment. 
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Appendix A 
Energy denominator transformations 

table iii 

Fundamental energy (frequency) difference contributions to homogeneously broadened (L) and inhomogeneously broadened 

(IB) electronic transitions up to second-order. 
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Appendix B 

Energy Denominators for Second-Order Processes 
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C. Sum and Difference Frequency Generation 
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